home *** CD-ROM | disk | FTP | other *** search
/ Libris Britannia 4 / science library(b).zip / science library(b) / UTILITIE / OKUMURA.ZIP / LZARI.C < prev    next >
Text File  |  1989-04-09  |  17KB  |  459 lines
  1. /**************************************************************
  2.         LZARI.C -- A Data Compression Program
  3.         (tab = 4 spaces)
  4. ***************************************************************
  5.         4/7/1989 Haruhiko Okumura
  6.         Use, distribute, and modify this program freely.
  7.         Please send me your improved versions.
  8.                 PC-VAN          SCIENCE
  9.                 NIFTY-Serve     PAF01022
  10.                 CompuServe      74050,1022
  11. **************************************************************/
  12. #include <stdio.h>
  13. #include <stdlib.h>
  14. #include <string.h>
  15. #include <ctype.h>
  16.  
  17. /********** Bit I/O **********/
  18.  
  19. FILE  *infile, *outfile;
  20. unsigned long int  textsize = 0, codesize = 0, printcount = 0;
  21.  
  22. void Error(char *message)
  23. {
  24.         printf("\n%s\n", message);
  25.         exit(EXIT_FAILURE);
  26. }
  27.  
  28. void PutBit(int bit)  /* Output one bit (bit = 0,1) */
  29. {
  30.         static unsigned int  buffer = 0, mask = 128;
  31.  
  32.         if (bit) buffer |= mask;
  33.         if ((mask >>= 1) == 0) {
  34.                 if (putc(buffer, outfile) == EOF) Error("Write Error");
  35.                 buffer = 0;  mask = 128;  codesize++;
  36.         }
  37. }
  38.  
  39. void FlushBitBuffer(void)  /* Send remaining bits */
  40. {
  41.         int  i;
  42.  
  43.         for (i = 0; i < 7; i++) PutBit(0);
  44. }
  45.  
  46. int GetBit(void)  /* Get one bit (0 or 1) */
  47. {
  48.         static unsigned int  buffer, mask = 0;
  49.  
  50.         if ((mask >>= 1) == 0) {
  51.                 buffer = getc(infile);  mask = 128;
  52.         }
  53.         return ((buffer & mask) != 0);
  54. }
  55.  
  56. /********** LZSS with multiple binary trees **********/
  57.  
  58. #define N                4096   /* size of ring buffer */
  59. #define F                  60   /* upper limit for match_length */
  60. #define THRESHOLD       2   /* encode string into position and length
  61.                                                    if match_length is greater than this */
  62. #define NIL                     N       /* index for root of binary search trees */
  63.  
  64. unsigned char  text_buf[N + F - 1];     /* ring buffer of size N,
  65.                         with extra F-1 bytes to facilitate string comparison */
  66. int             match_position, match_length,  /* of longest match.  These are
  67.                         set by the InsertNode() procedure. */
  68.                 lson[N + 1], rson[N + 257], dad[N + 1];  /* left & right children &
  69.                         parents -- These constitute binary search trees. */
  70.  
  71. void InitTree(void)  /* Initialize trees */
  72. {
  73.         int  i;
  74.  
  75.         /* For i = 0 to N - 1, rson[i] and lson[i] will be the right and
  76.            left children of node i.  These nodes need not be initialized.
  77.            Also, dad[i] is the parent of node i.  These are initialized to
  78.            NIL (= N), which stands for 'not used.'
  79.            For i = 0 to 255, rson[N + i + 1] is the root of the tree
  80.            for strings that begin with character i.  These are initialized
  81.            to NIL.  Note there are 256 trees. */
  82.  
  83.         for (i = N + 1; i <= N + 256; i++) rson[i] = NIL;       /* root */
  84.         for (i = 0; i < N; i++) dad[i] = NIL;   /* node */
  85. }
  86.  
  87. void InsertNode(int r)
  88.         /* Inserts string of length F, text_buf[r..r+F-1], into one of the
  89.            trees (text_buf[r]'th tree) and returns the longest-match position
  90.            and length via the global variables match_position and match_length.
  91.            If match_length = F, then removes the old node in favor of the new
  92.            one, because the old one will be deleted sooner.
  93.            Note r plays double role, as tree node and position in buffer. */
  94. {
  95.         int  i, p, cmp, temp;
  96.         unsigned char  *key;
  97.  
  98.         cmp = 1;  key = &text_buf[r];  p = N + 1 + key[0];
  99.         rson[r] = lson[r] = NIL;  match_length = 0;
  100.         for ( ; ; ) {
  101.                 if (cmp >= 0) {
  102.                         if (rson[p] != NIL) p = rson[p];
  103.                         else {  rson[p] = r;  dad[r] = p;  return;  }
  104.                 } else {
  105.                         if (lson[p] != NIL) p = lson[p];
  106.                         else {  lson[p] = r;  dad[r] = p;  return;  }
  107.                 }
  108.                 for (i = 1; i < F; i++)
  109.                         if ((cmp = key[i] - text_buf[p + i]) != 0)  break;
  110.                 if (i > THRESHOLD) {
  111.                         if (i > match_length) {
  112.                                 match_position = (r - p) & (N - 1);
  113.                                 if ((match_length = i) >= F) break;
  114.                         } else if (i == match_length) {
  115.                                 if ((temp = (r - p) & (N - 1)) < match_position)
  116.                                         match_position = temp;
  117.                         }
  118.                 }
  119.         }
  120.         dad[r] = dad[p];  lson[r] = lson[p];  rson[r] = rson[p];
  121.         dad[lson[p]] = r;  dad[rson[p]] = r;
  122.         if (rson[dad[p]] == p) rson[dad[p]] = r;
  123.         else                   lson[dad[p]] = r;
  124.         dad[p] = NIL;  /* remove p */
  125. }
  126.  
  127. void DeleteNode(int p)  /* Delete node p from tree */
  128. {
  129.         int  q;
  130.  
  131.         if (dad[p] == NIL) return;  /* not in tree */
  132.         if (rson[p] == NIL) q = lson[p];
  133.         else if (lson[p] == NIL) q = rson[p];
  134.         else {
  135.                 q = lson[p];
  136.                 if (rson[q] != NIL) {
  137.                         do {  q = rson[q];  } while (rson[q] != NIL);
  138.                         rson[dad[q]] = lson[q];  dad[lson[q]] = dad[q];
  139.                         lson[q] = lson[p];  dad[lson[p]] = q;
  140.                 }
  141.                 rson[q] = rson[p];  dad[rson[p]] = q;
  142.         }
  143.         dad[q] = dad[p];
  144.         if (rson[dad[p]] == p) rson[dad[p]] = q;
  145.         else                   lson[dad[p]] = q;
  146.         dad[p] = NIL;
  147. }
  148.  
  149. /********** Arithmetic Compression **********/
  150.  
  151. /*  If you are not familiar with arithmetic compression, you should read
  152.                 I. E. Witten, R. M. Neal, and J. G. Cleary,
  153.                         Communications of the ACM, Vol. 30, pp. 520-540 (1987),
  154.         from which much have been borrowed.  */
  155.  
  156. #define M   15
  157.  
  158. /*      Q1 (= 2 to the M) must be sufficiently large, but not so
  159.         large as the unsigned long 4 * Q1 * (Q1 - 1) overflows.  */
  160.  
  161. #define Q1  (1UL << M)
  162. #define Q2  (2 * Q1)
  163. #define Q3  (3 * Q1)
  164. #define Q4  (4 * Q1)
  165. #define MAX_CUM (Q1 - 1)
  166.  
  167. #define N_CHAR  (256 - THRESHOLD + F)
  168.         /* character code = 0, 1, ..., N_CHAR - 1 */
  169.  
  170. unsigned long int  low = 0, high = Q4, value = 0;
  171. int  shifts = 0;  /* counts for magnifying low and high around Q2 */
  172. int  char_to_sym[N_CHAR], sym_to_char[N_CHAR + 1];
  173. unsigned int
  174.         sym_freq[N_CHAR + 1],  /* frequency for symbols */
  175.         sym_cum[N_CHAR + 1],   /* cumulative freq for symbols */
  176.         position_cum[N + 1];   /* cumulative freq for positions */
  177.  
  178. void StartModel(void)  /* Initialize model */
  179. {
  180.         int ch, sym, i;
  181.  
  182.         sym_cum[N_CHAR] = 0;
  183.         for (sym = N_CHAR; sym >= 1; sym--) {
  184.                 ch = sym - 1;
  185.                 char_to_sym[ch] = sym;  sym_to_char[sym] = ch;
  186.                 sym_freq[sym] = 1;
  187.                 sym_cum[sym - 1] = sym_cum[sym] + sym_freq[sym];
  188.         }
  189.         sym_freq[0] = 0;  /* sentinel (!= sym_freq[1]) */
  190.         position_cum[N] = 0;
  191.         for (i = N; i >= 1; i--)
  192.                 position_cum[i - 1] = position_cum[i] + 10000 / (i + 200);
  193.                         /* empirical distribution function (quite tentative) */
  194.                         /* Please devise a better mechanism! */
  195. }
  196.  
  197. void UpdateModel(int sym)
  198. {
  199.         int i, c, ch_i, ch_sym;
  200.  
  201.         if (sym_cum[0] >= MAX_CUM) {
  202.                 c = 0;
  203.                 for (i = N_CHAR; i > 0; i--) {
  204.                         sym_cum[i] = c;
  205.                         c += (sym_freq[i] = (sym_freq[i] + 1) >> 1);
  206.                 }
  207.                 sym_cum[0] = c;
  208.         }
  209.         for (i = sym; sym_freq[i] == sym_freq[i - 1]; i--) ;
  210.         if (i < sym) {
  211.                 ch_i = sym_to_char[i];    ch_sym = sym_to_char[sym];
  212.                 sym_to_char[i] = ch_sym;  sym_to_char[sym] = ch_i;
  213.                 char_to_sym[ch_i] = sym;  char_to_sym[ch_sym] = i;
  214.         }
  215.         sym_freq[i]++;
  216.         while (--i >= 0) sym_cum[i]++;
  217. }
  218.  
  219. static void Output(int bit)  /* Output 1 bit, followed by its complements */
  220. {
  221.         PutBit(bit);
  222.         for ( ; shifts > 0; shifts--) PutBit(! bit);
  223. }
  224.  
  225. void EncodeChar(int ch)
  226. {
  227.         int  sym;
  228.         unsigned long int  range;
  229.  
  230.         sym = char_to_sym[ch];
  231.         range = high - low;
  232.         high = low + (range * sym_cum[sym - 1]) / sym_cum[0];
  233.         low +=       (range * sym_cum[sym    ]) / sym_cum[0];
  234.         for ( ; ; ) {
  235.                 if (high <= Q2) Output(0);
  236.                 else if (low >= Q2) {
  237.                         Output(1);  low -= Q2;  high -= Q2;
  238.                 } else if (low >= Q1 && high <= Q3) {
  239.                         shifts++;  low -= Q1;  high -= Q1;
  240.                 } else break;
  241.                 low += low;  high += high;
  242.         }
  243.         UpdateModel(sym);
  244. }
  245.  
  246. void EncodePosition(int position)
  247. {
  248.         unsigned long int  range;
  249.  
  250.         range = high - low;
  251.         high = low + (range * position_cum[position    ]) / position_cum[0];
  252.         low +=       (range * position_cum[position + 1]) / position_cum[0];
  253.         for ( ; ; ) {
  254.                 if (high <= Q2) Output(0);
  255.                 else if (low >= Q2) {
  256.                         Output(1);  low -= Q2;  high -= Q2;
  257.                 } else if (low >= Q1 && high <= Q3) {
  258.                         shifts++;  low -= Q1;  high -= Q1;
  259.                 } else break;
  260.                 low += low;  high += high;
  261.         }
  262. }
  263.  
  264. void EncodeEnd(void)
  265. {
  266.         shifts++;
  267.         if (low < Q1) Output(0);  else Output(1);
  268.         FlushBitBuffer();  /* flush bits remaining in buffer */
  269. }
  270.  
  271. int BinarySearchSym(unsigned int x)
  272.         /* 1      if x >= sym_cum[1],
  273.            N_CHAR if sym_cum[N_CHAR] > x,
  274.            i such that sym_cum[i - 1] > x >= sym_cum[i] otherwise */
  275. {
  276.         int i, j, k;
  277.  
  278.         i = 1;  j = N_CHAR;
  279.         while (i < j) {
  280.                 k = (i + j) / 2;
  281.                 if (sym_cum[k] > x) i = k + 1;  else j = k;
  282.         }
  283.         return i;
  284. }
  285.  
  286. int BinarySearchPos(unsigned int x)
  287.         /* 0 if x >= position_cum[1],
  288.            N - 1 if position_cum[N] > x,
  289.            i such that position_cum[i] > x >= position_cum[i + 1] otherwise */
  290. {
  291.         int i, j, k;
  292.  
  293.         i = 1;  j = N;
  294.         while (i < j) {
  295.                 k = (i + j) / 2;
  296.                 if (position_cum[k] > x) i = k + 1;  else j = k;
  297.         }
  298.         return i - 1;
  299. }
  300.  
  301. void StartDecode(void)
  302. {
  303.         int i;
  304.  
  305.         for (i = 0; i < M + 2; i++)
  306.                 value = 2 * value + GetBit();
  307. }
  308.  
  309. int DecodeChar(void)
  310. {
  311.         int      sym, ch;
  312.         unsigned long int  range;
  313.  
  314.         range = high - low;
  315.         sym = BinarySearchSym((unsigned int)
  316.                 (((value - low + 1) * sym_cum[0] - 1) / range));
  317.         high = low + (range * sym_cum[sym - 1]) / sym_cum[0];
  318.         low +=       (range * sym_cum[sym    ]) / sym_cum[0];
  319.         for ( ; ; ) {
  320.                 if (low >= Q2) {
  321.                         value -= Q2;  low -= Q2;  high -= Q2;
  322.                 } else if (low >= Q1 && high <= Q3) {
  323.                         value -= Q1;  low -= Q1;  high -= Q1;
  324.                 } else if (high > Q2) break;
  325.                 low += low;  high += high;
  326.                 value = 2 * value + GetBit();
  327.         }
  328.         ch = sym_to_char[sym];
  329.         UpdateModel(sym);
  330.         return ch;
  331. }
  332.  
  333. int DecodePosition(void)
  334. {
  335.         int position;
  336.         unsigned long int  range;
  337.  
  338.         range = high - low;
  339.         position = BinarySearchPos((unsigned int)
  340.                 (((value - low + 1) * position_cum[0] - 1) / range));
  341.         high = low + (range * position_cum[position    ]) / position_cum[0];
  342.         low +=       (range * position_cum[position + 1]) / position_cum[0];
  343.         for ( ; ; ) {
  344.                 if (low >= Q2) {
  345.                         value -= Q2;  low -= Q2;  high -= Q2;
  346.                 } else if (low >= Q1 && high <= Q3) {
  347.                         value -= Q1;  low -= Q1;  high -= Q1;
  348.                 } else if (high > Q2) break;
  349.                 low += low;  high += high;
  350.                 value = 2 * value + GetBit();
  351.         }
  352.         return position;
  353. }
  354.  
  355. /********** Encode and Decode **********/
  356.  
  357. void Encode(void)
  358. {
  359.         int  i, c, len, r, s, last_match_length;
  360.  
  361.         fseek(infile, 0L, SEEK_END);
  362.         textsize = ftell(infile);
  363.         if (fwrite(&textsize, sizeof textsize, 1, outfile) < 1)
  364.                 Error("Write Error");  /* output size of text */
  365.         codesize += sizeof textsize;
  366.         if (textsize == 0) return;
  367.         rewind(infile);  textsize = 0;
  368.         StartModel();  InitTree();
  369.         s = 0;  r = N - F;
  370.         for (i = s; i < r; i++) text_buf[i] = ' ';
  371.         for (len = 0; len < F && (c = getc(infile)) != EOF; len++)
  372.                 text_buf[r + len] = c;
  373.         textsize = len;
  374.         for (i = 1; i <= F; i++) InsertNode(r - i);
  375.         InsertNode(r);
  376.         do {
  377.                 if (match_length > len) match_length = len;
  378.                 if (match_length <= THRESHOLD) {
  379.                         match_length = 1;  EncodeChar(text_buf[r]);
  380.                 } else {
  381.                         EncodeChar(255 - THRESHOLD + match_length);
  382.                         EncodePosition(match_position - 1);
  383.                 }
  384.                 last_match_length = match_length;
  385.                 for (i = 0; i < last_match_length &&
  386.                                 (c = getc(infile)) != EOF; i++) {
  387.                         DeleteNode(s);  text_buf[s] = c;
  388.                         if (s < F - 1) text_buf[s + N] = c;
  389.                         s = (s + 1) & (N - 1);
  390.                         r = (r + 1) & (N - 1);
  391.                         InsertNode(r);
  392.                 }
  393.                 if ((textsize += i) > printcount) {
  394.                         printf("%12ld\r", textsize);  printcount += 1024;
  395.                 }
  396.                 while (i++ < last_match_length) {
  397.                         DeleteNode(s);
  398.                         s = (s + 1) & (N - 1);
  399.                         r = (r + 1) & (N - 1);
  400.                         if (--len) InsertNode(r);
  401.                 }
  402.         } while (len > 0);
  403.         EncodeEnd();
  404.         printf("In : %lu bytes\n", textsize);
  405.         printf("Out: %lu bytes\n", codesize);
  406.         printf("Out/In: %.3f\n", (double)codesize / textsize);
  407. }
  408.  
  409. void Decode(void)
  410. {
  411.         int  i, j, k, r, c;
  412.         unsigned long int  count;
  413.  
  414.         if (fread(&textsize, sizeof textsize, 1, infile) < 1)
  415.                 Error("Read Error");  /* read size of text */
  416.         if (textsize == 0) return;
  417.         StartDecode();  StartModel();
  418.         for (i = 0; i < N - F; i++) text_buf[i] = ' ';
  419.         r = N - F;
  420.         for (count = 0; count < textsize; ) {
  421.                 c = DecodeChar();
  422.                 if (c < 256) {
  423.                         putc(c, outfile);  text_buf[r++] = c;
  424.                         r &= (N - 1);  count++;
  425.                 } else {
  426.                         i = (r - DecodePosition() - 1) & (N - 1);
  427.                         j = c - 255 + THRESHOLD;
  428.                         for (k = 0; k < j; k++) {
  429.                                 c = text_buf[(i + k) & (N - 1)];
  430.                                 putc(c, outfile);  text_buf[r++] = c;
  431.                                 r &= (N - 1);  count++;
  432.                         }
  433.                 }
  434.                 if (count > printcount) {
  435.                         printf("%12lu\r", count);  printcount += 1024;
  436.                 }
  437.         }
  438.         printf("%12lu\n", count);
  439. }
  440.  
  441. int main(int argc, char *argv[])
  442. {
  443.         char  *s;
  444.  
  445.         if (argc != 4) {
  446.                 printf("'lzari e file1 file2' encodes file1 into file2.\n"
  447.                            "'lzari d file2 file1' decodes file2 into file1.\n");
  448.                 return EXIT_FAILURE;
  449.         }
  450.         if ((s = argv[1], s[1] || strpbrk(s, "DEde") == NULL)
  451.          || (s = argv[2], (infile  = fopen(s, "rb")) == NULL)
  452.          || (s = argv[3], (outfile = fopen(s, "wb")) == NULL)) {
  453.                 printf("??? %s\n", s);  return EXIT_FAILURE;
  454.         }
  455.         if (toupper(*argv[1]) == 'E') Encode();  else Decode();
  456.         fclose(infile);  fclose(outfile);
  457.         return EXIT_SUCCESS;
  458. }
  459.